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Abstract 

A second-order form of discrete Kalman filtering equations is proposed as a candidate state estimator for 
efficient simulations of control-structure interactions in coupled physical coordinate configurations as opposed 
to decoupled modal coordinates. The resulting matrix equation of the present state estimator consists of the 
same symmetric, sparse N x N coupled matrices of the governing structural dynamics equations as opposed 
to unsymmetric 2N x 2 N state space-based estimators. Thus, in addition to substantial computational 
efficiency improvement, the present estimator can be applied to control-structure design optimization for 
which the physical coordinates associated with the mass, damping and stiffness matrices of the structure are 
needed instead of modal coordinates. 
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Introduction 


Current practice in the design, modeling and analysis of flexible large space structures is by and large 
based on the finite element method and the associated software. The resulting discrete equations of motion 
for structures, both in terms of physical coordinates and of modal coordinates, are expressed in a second- 
order form. As a result, the structural engineering community has been investing a considerable amount of 
research and development resources to develop computer-oriented discrete modeling tools, analysis methods 
and interface capabilities with design synthesis procedures; all of these exploiting the characteristics of 
second-order models. Recent work in the area of structural dynamics simulation and massively-parallel 
processing also rely on the second-order equation forms. 

On the other hand, modern linear control theory has its roots firmly in a first-order form of the governing 
differential equations, e.g., Kwakernaak and Sivan 1 . Thus, several investigators have addressed the issues of 
interfacing second-order structural systems and control theory based on the first-order form . As a result 
of these studies, it has become straightforward for one to synthesize direct state feedback based control laws 
within the framework of a first-order control theory and then to recast the resulting control laws in terms of 
the second-order structural systems. 

Unfortunately, controllers based on a first-order state estimator are difficult to express in a pure second-order 
form because the first-order estimator implicitly incorporates an additional filter equation 7 . However a recent 
work by Juang and Maghami 8 has enabled the first-order filter gain matrices to be synthesized using only 
second-order equations. To complement the second-order gain synthesis, the objective of the present paper 
is to develop a second-order based simulation procedure for first-order estimators. The particular class of 
first-order dynamic compensation chosen for study are the Kalman Filter based state estimators as applied to 
second-order structural systems. The proposed procedure permits simulation of first-order estimators with 
nearly the same solution procedure used for treating the structural dynamics equation. Hence, the reduced 
size of system matrices and the computational techniques that are tailored to sparse second-order structural 
systems may be employed. As will be shown, the proposed procedure hinges on discrete time integration 
formulas to effectively reduce the continuous time Kalman Filter to a set of second-order difference equations. 

The primary goal of the proposed procedure is the incorporation of this general form of state estimation 
as a simulation tool in partitioned control-structure interaction (CSI) analyses. It is expected that Kalman 
Filters for real-time control of linear time-invariant systems would be implemented in the most efficient form 
available, typically a real mode-decoupled state space realization. For analytical studies of CSI systems, 
however, where the objective is frequently simultaneous optimization of controls and structures as in Belvin , 
the use of such a modal form must be weighed against the preprocessing tasks required to generate the model. 
In these cases, much more flexible controllers expressed in terms of the physical coordinates instead of the 
modal coordinates are sought, ones which can readily adapt to iterative changes in the structural parameters. 
One such control law synthesis has been proposed and demonstrated to be effective for CSI optimization . 
These studies did not account, however, for dynamic compensation when full state feedback control was 
utilized. With the discrete Kalman Filter proposed herein, a general form of dynamic compensation can be 
integrated into CSI simulation and optimization which does not impose limits on the designs of the feedback 
gains or the filter gains. 

The paper first reviews of the conventional first-order representation of the continuous second-order structural 
equations of motion, in which the state variables are defined as the displacements variables i of the second- 
order structural model and the velocities i. An examination of the corresponding first-order Kalman filtering 
equations indicates that, due to the difference in the derivative of the estimated displacement (£i) and the 
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estimated velocity (i), transformation of the first-order estimator into an equivalent second-order estimator 
requires the time derivative of measurement data, a process not recommended for practical implementation. 

Next, a transformation via a generalized momentum is introduced to recast the structural equations of motion 
in a general first-order setting. It is shown that discrete time numerical integration followed by reduction 
of the resulting difference equations circumvents the need for the time derivative of measurements to solve 
Kalman filtering equations in a second-order framework. Hence, the Kalman filter equations can be solved 
using a second-order solution software package. 

Subsequently, computer implementation aspects of the present second-order estimator are presented. Several 
computational paths are discussed in the context of discrete and continuous time simulation. For continuous 
time control simulation, an equation augmentation is introduced to exploit the symmetry and sparsity of 
the attendant matrices by maintaining state dependent control and observer terms on the right-hand-side 
(RHS) of the filter equations. In addition, the computational efficiency of the present second-order filter as 
compared to the first-order form is presented. 


Continuous Formulation of State 
Estimators for Structural Systems 

Linear, second-order discrete structural models can be expressed as 

Mi 4 Dx -f Kx — Bu + Gw , x(0) = xo , i(0) = io (!) 

U = —Z\T — Z 2 X 

with the associated measurements 


z = H X X 4 H 2 x 4 v 


( 2 ) 


where M, D , K are the mass damping and stiffness matrices of size ( N x N)] x is the structural displacement 
vector, {N x 1); u is the active control force (m x 1); B is a constant force distribution matrix ( N x m); z 
is a set of measurements (r x 1); H x and H 2 are the measurement distribution matrices (r x N)\ Z x and 
Z 2 are the control feedback gain matrices (m x N)\ w and v are zero-mean, white Gaussian processes with 
their respective covariances Q and R\ and the superscript dot designates time differentiation. In the present 
study, we will restrict ourselves to the case wherein Q and R are uncorrelated with each other and the initial 
conditions x 0 and io are also themselves jointly Gaussian with known means and covariances. 

The conventional representation of (1) in a first-order form is facilitated by 

(x x -x 

J x 2 — x — X\ (3) 

l Mx 2 — Mi = Bu + Gw — Dx 2 — Kx 1 


which, when cast in a first-order form, can be expressed as 

( Eq — Fq 4 Bu 4 Gw, q = ( x x x 2 ) T 
\z = Hq-t v 


where 
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( 5 ) 


B = 




It is well-known that the Kalman filtering equations 10 ' 11 for (4) can be shown to be (see, e.g., Arnold and 
Laub 3 ): 


where 


z 


Eq = Fq + Bu + EPH T R 


z — Hq , 



L ' 


9 = 


- l z 

( 6 ) 

( Xy 1 f X 1 

\ X 2 / ~~ \ x f 

( 7 ) 


in which U and L are positive definite matrices, q is the state estimation vector, and the matrix P is 
determined by the Riccati equation 1,3 

EPE t = FPE t -f EPF t - EPH T R~ l HPE T + GQG t (8) 


The inherent difficulty of reducing the first-order Kalman filtering equations given by (6) to second order 
form can be appreciated if one attempts to write (6) in a form introduced in (3): 

{ a) zi=x 

b) x 2 — x — i\ — Lyz (9) 

c) Mx 2 = — Dx 2 — Kx\ + Bu + M L 2 z 

where 

L x = ( HyU + H 2 S) T R~\ L 2 = ( HyS T + H 2 L) T IT 1 

Note from (9b) that i 2 / i\. In other words, the time derivative of the estimated displacement (£) is not 
the same as the estimated velocity (i); hence, xy and x 2 must be treated as two independent variables, an 
important observation somehow overlooked in Hashemipour and Laub 12 . 

Of course, although not practical, one can eliminate x 2 from (9). Assuming and x 2 are differentiable, 
differentiate (9b) and multiply both sides by M to obtain 

Mi\ = Mx 2 + MLyz (10) 

Substituting Mx 2 from (9c) and x 2 from (9b) in (10) yields 

Mxy — —D(xy — L\z) — Kxy + Bu + M L 2 Z + M LyZ (11) 


which, upon rearrangements, becomes 

Mxy + Dx 1 4* Kx 1 = Bu + M L 2 z + MLyz -f DLyz (12) 

There are two difficulties with the above second-order estimator. First, the numerical solution of (12) 
involves the computation of when rate measurements are made. The accuracy of this computation is in 
general very susceptible to errors caused in numerical differentiation of xy. Second, and most important, the 
numerical evaluation of z that is required in (12) assumes that the derivative of measurement information is 
available which should be avoided in practice. We now present a computational procedure that circumvents 
the need for computing measurement derivatives and that enables one to construct estimators based on the 
second-order model form. 
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Second-Order Transformation of 
Continuous Kalman Filtering Equations 

This section presents a transformation of the continuous time first-order Kalman filter to a discrete time set of 
second-order difference equations for digital implementation. The procedure avoids the need for measurement 
derivative information. In addition, the sparsity and symmetry of the original mass, damping and stiffness 
matrices can be maintained. Prior to describing the numerical integration procedure, a transformation based 
on generalized momenta is presented which is later used to improve computational efficiency of the equation 
solution. 


Generalized Momenta 

Instead of the conventional transformation (3) of the second-order structural system (1) into a first-order 
form, let us consider the following generalized momenta (see, e.g., Jensen 13 and Felippa and Park 14 ): 

( a) Xl =x ( 13 ) 

\ b) x 2 = AMix+Cxt 

where A and C are constant matrices to be chosen. Note that AM should be nonsingular in order to obtain 
an equivalent form of (1). Time differentiation of (13b) yields 

i2 = AM xy + Cx\ (14) 


Substituting (1) via (13a) into (14), one obtains 

X 2 — A(Bu + Gw) — (AD — C)i\ — AKx\ 

Finally, pairing of (13b) and (15) gives the following first-order form: 

. r °UM + [ c - J 1M = 

[AD-C I \ \i 2 J + [AK 0 J \*2 J 

[ 0 1 

[A(Bu + Gw) \ 

The associated Kalman filtering equation can be shown to be of the following form: 

[ AM °UM+f c - 7 U fi W 0 U 

AD-C /J \*2 J [AK 0 J 1 *2 J \ ABu j 

[ AM 0] [Li] _ 

[AD-C I J [l 2 J Z 

where 

Li = (HiU + H 2 S) T R~ l , L 2 = (HiS t + H 2 L) T R- 1 

and Hi and H 2 correspond to a modified form of measurements expressed as 


z = H\X + H 2 x = H 1X1 + H 2 x i 


(15) 


(16) 


(17) 


( 18 ) 


where 

Hi = Hi - H 2 M~ l A~ l C, n 2 = H 2 M~ l A - 1 
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Clearly, as in the conventional first-order form (9), fi and x-i in (17) are now two independent variables. 
Specifically, the case of A — M~ l and C = 0 corresponds to (3) with *2 = *1- However, as we shall see 
below, the Kalman filtering equations based on the generalized momenta (13) offer several computational 
advantages over (3). 


Numerical Integration 

At this juncture it is noted that in the previous section one first performs the elimination of ij in order 
to obtain a second-order equation, then performs the numerical solution of the resulting equation. This 
approach has the disadvantage of having to deal with the time derivative of measurement data. To avoid 
this, we will first integrate numerically the associated Kalman filtering equation (17). 

The direct time integration formula we propose to employ is a mid-point version of the trapezoidal rule: 





1 n + l/2 

[I j +6 {l 

} 

fi 1 r +1/2 f 

1 ) 

ill" 

Ul -{ 



(19) 


where the superscript n denotes the discrete time interval t n = nh , h is the time increment and 6 — h/2. 

It should be noted that we have chosen the trapezoidal rule due to its unconditional stability and high 
accuracy while it does not introduce any numerical damping (see, for example, Dahlquist 15 and Park 16 ). 
Contamination of damping from numerical dissipation can not only adversely affect the solution accuracy 
but lead to misinterpretation of the simulation results. 


Time discretization of (17) by (19a) at the n -f 1/2 time step yields 


AM 

AD-C 


01 fi" +1/2 -z? 

/ J U2 + 1/2 -x5 




AM 0 
AD-C I 



f " + 1/2 + 6 { ABu n+1 ' 2 } 


( 20 ) 


The above difference equations require the solution of matrix equations of 2N variables, namely, in terms 
of the two variables x^ lf2 and x7 +1/2 , each with a size of N. To reduce the above coupled equations of 
order 2 N into the corresponding ones of order N y we proceed in the following way by exploiting the nature 
of parametric matrices of A and C as introduced in (13). To this end, we write out (20) as two coupled 
difference equations as follows: 

AM(x n ^ 1 ' 2 - i?) + «(Ci" +1/2 - x" +1/2 ) 

= 6AML l z n+l/2 (21) 

(AD - C)(x^ +1/2 - x?) + (x^ +1/2 - x$) + 6AKx n + 1 ' 2 

= S(AD - C)L,f n+1/2 + 6L 7 z n+1/2 + 6ABu n+l > 2 (22) 

Multiplying (22) by S and adding the resulting equation to (21) yields 

A(M + 6D + 6 2 /OxJ +1/2 = (AM + 5(AD - C ))x? + 6x” 

+{6AMli + 6 2 (AD - C)L X + 6 2 L 7 }z n+1/ 2 + 6 2 ABu n+l/ 2 (23) 
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Of several possible choices for matrices A and 5, we will examine the following two specific cases: 


f a) A — I t C=D 

{ b) A — M~ l , C = 0 


( 24 ) 


where the mass matrix M is nonsingular due to its physically positive definite nature since the kinetic energy 
of structural system is positive for any admissable motion. It is noted that the above two choices, although 
mathematically equivalent, lead to different computational implementations as discussed below. 

The choice of (24a) reduces (23) to: 


(M + 6D + 6 2 K)x n ^ 1 ' 2 = Mx'l + 6% + 6 2 Bu n+1/ 2 
+6{ML l + 6L 2 )z n+l/2 

so that once x” +1 ^ 2 is computed, x, + 1 ^ 2 is obtained from (22) rewritten as 

x£ +1/2 = ±1 + 6g n - 6Kx n + 1 ' 2 


(25) 


(26) 


u/nArp 

g n = Bu n+I/2 + L 2 £ n+1/2 


(27) 


which is already computed in order to construct the right-hand side of (25). Hence, Kx ” +1 ^ 2 is the only 
additional computation needed to obtain x ^ 1 ^ 2 . It is noted that neither any numerical differentiation nor 
matrix inversion is required in computing x!^ 1 ^ 2 . This has been achieved through the introduction of the 
general transformation (13) and the particular choice of the parameter matrices given by (24a). 

On the other hand, if one chooses the conventional representation (24b), the solution of x " +1 ^ 2 is obtained 
from (23) 

(M + 6D + 6 2 K )x " +1/2 = {M + 6D)x ? + 6Mx £ 

+6{(M + SD)^ 4- SMU}^ 2 + * 2 Bu n+l/2 (28) 

Once x” +1/r2 is obtained, i^ 1 ^ 2 can be computed either by 


-n+l/2 

X 2 


— V X 1 


■*?)/*- L,z 


-n + l/2 


(29) 


;n+l/2 


which is not accurate due to the numerical differentiation to obtain x t , or by (22) 


X 2 +1/2 = xl + 6g n - 6M~ x Kx[ +l ' 2 - 

M~ 1 D(x" +1/2 - x?) + 6M~ l DL l z n+1/2 (30) 


which involves two additional matrix-vector multiplications, when D ^ 0, as compared with the choice of 
A = I and C — D. Thus (24a) is the preferred representation in a first-order form of the second-order 
structural dynamics equations (1) and is used in the remainder of this work. 


Decoupling Of Difference Equations 

We have seen in the previous section, instead of solving the first-order Kalman filtering equations of 2n 
variables for the structural dynamics systems (1), the solution of the implicit time-discrete estimator equation 
(25) of n variables can potentially offer a substantial computational saving by exploiting the reduced size 
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and sparsity of M , D and I < . This assumes that z n+1/2 and u n+1/2 are available, which is not the case since 
at the n th time step 


u "+l/2 


Z n + 1/2 _ z 


= -Z l x" +in *- -Z 2 x 2 +l/2 

(31) 

n+l ' 2 - H ,i? +l/2 - H 2 x n 2 +l ' 2 

(32) 


requires both x” +1/2 and x!J +1/2 even if z n+1/2 is assumed to be known from measurements or by solution 
of (1). Note in (32), the control gain matrices are transformed by 

Zi = Z x - Z 2 M- l A~ l C } Z 2 = Z 2 M- l A- 1 


There are two distinct approaches to decouple (25) and (26) as described in the following sections. 


Discrete Time Update 

For systems utilizing discrete-time (i.e. sample and hold) control, equations (31) and (32) become 

w n+l/2 ~ _ z 2 x^ (33) 

f n+1 / 2 ~ z n - H x x n Y - H 2 X 2 (34) 

The time integration step size of the estimator must then be equal to the sample rate of the control, while the 
continuous structural equations may also be integrated at the same rate or at some fraction of the sampling 
rate for simulation accuracy considerations. For the present purposes, we have assumed that the sampling 
interval is the same as the integration time stepsize. 

Discrete time simulation is quite simple to implement as the control force and state corrections are treated 
with no approximation on the right-hand-side (RHS) of (25) and (26). Should continuous time simulation 
be required, a different approach is necessary. 


Continuous Time Update 


To simulate the system given in (25) and (26) in continuous time, strictly speaking, one must rearrange 
(25) and (26) so that the terms involving x” +1/2 and x 2 + i/2 are augmented to the left-hand-side (LHS) of 
the equations. However, this augmentation into the solution matrix (M + 6D + S 2 K ) would destroy the 
computational advantages of the matrix sparsity and symmetry. Thus, a partitioned solution procedure has 
been developed for continuous time simulation as described in Park and Belvin 17 . The procedure, briefly 
outlined herein, maintains the control force and state correction on the RHS of the equations as follows. 

First, i" +1 ^ 2 and x j +1 ^ 2 are predicted by 


-n+l/2 

x lp 


-1 • 


-n+l/2 

x 2p 


(35) 


However, instead of direct substitution of the above predicted quantity to obtain tip ^ and ip based 

, , . . f n + l/2 , -n + l/2 

on (31) and (32), equation augmentations are introduced to improve the accuracy ot u p ana z p 
Of several augmentation procedures that are applicable to construct discrete filters for the computations of 
u”+i/2 and 2 n+1 / 2 , we substitute (26) into (31) and (32) to obtain 


{ u n+l/2 _ ^2i*i +l/9 - Z 2 {il - 6Kx" +l,2 + 

6Bun+'/ 2 + 6Z, 2 f n+1 /2) 

£" + 1/2 = *" + 1/2 _ fliij +1/2 — 

H 2 {x 5 - 6Ki" +U2 + 6Bu n+1 ' 2 + 6L 2 z n+l/2 ) 
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Rearranging the above coupled equations, one obtains 


'(I + 6Z 2 B) 6Z 2 L 2 lf«" +1/2 \ = 

6H 2 B (I + 6H 2 L 2 )\ \z" +,/2 J 

/ -Z 2 x n 2 -(Zi-f>Z 2 K)x \^ 12 \ (m 

\ z "+i/2 _ h 2 x J - (Hi - SH 2 K)xi +l/2 J 

which corresponds to a first order filter to reduce the errors in computing x 2 = Mx + Dx. A second-order 
discrete filter for computing u and z can be obtained by differentiating u and z to obtain 

f u = —Z\i\ - Z 2 x 2 (38) 

\z = z — H\i\ - H 2 X 2 


and then substituting and x 2 from (17). Subsequently, (19) is applied to integrate the equations for u 
and z which yields 


7 + SZ 2 B + S 7 Z\M~ l B 6(Z 2 L 2 + Zih + 6Z t M- l L 
6(H 2 B + 6H l M~ l B) I + 6H l {L x +6M- iL L 2 ) + SH 

f 11 " 1 r Z i M- l (xq - 6ivx” +1/2 - Dx n ^ l/2 ) + Z 2 Kx n l + l 


(39) 


The net effects of this augmentation are to filter out the errors committed in estimating both x x and x 2 . 
Solution of (39) for u n+1 ^ 2 and z n + 1 / 2 permits (25) and (26) to be solved in continuous time for x" * 1/2 and 
^n+i /2 g u b S equently, (19b) is used for x"* 1 and xj* 1 . 


The preceding augmentation (39) leads to an accurate estimate of the control force and state estimation error 
correction at the (n-t-1/2) time step. Although (39) involves the solution of an additional algebraic equation, 
the equation size is relatively small ( size = number of actuators (m) plus the number of measurements (r) 
). Thus, (39) is an efficient method for continuous time simulation of the Kalman filter equations provided 
the size of (39) is significantly lower than the first order form of (4). The next section discusses the relative 
efficiency of the present method and the conventional first order solution. More details on the equation 
augmentation procedure (39) may be found in Park and Belvin 17 . 

Finally, it is noted that by following a similar time discretization procedure adopted for computing x l 
and x" +l/ \ the structural dynamics equation (1) can be solved by 


/ (Af + 6D + 6 2 /C)x” +1/2 = Afx” + 6x” + 6 7 Bu n+l f 7 (AQ) 

\ X 2 +1/2 = *3 + 6Bu n + l f 2 - 6Kx" +l/2 

Thus, numerical solutions of the structural dynamics equation (1) and the filter equation (20) can be carried 
out within the second-order solution context, thus realizing substantial computational simplicity compared 
with the solution of first-order systems of equations (4) and the corresponding first-order estimator equations 
( 6 ). 

It is emphasized that the solutions of both the structural displacement x and the reconstructed displacement 
x employ the same solution matrix, (M + 6D 4* <5 2 /f). The computational stability of the present procedure 
can be examined as investigated in Park 18 and Park and Felippa 19,20 . The result, when applied to the present 
case, can be stated as 

6 2 Amax < 1 ( 41 ) 
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where A ma jc is the maximum eigenvalue of 

(A 2 1 + A Z 2 B + Z\M~ l B)y = 0 (42) 

Typically the control laws are formulated in terms of low-frequency response components, viz., 

BcxG T KG (43) 

for the displacement feedback case where G is a projection matrix that extracts only low-frequency compo- 
nents from the structural stiffness matrix. Hence, A max is in practice several orders of magnitude smaller 
than fimax of the structural dynamics eigenvalue problem: 

\i My-Ky (44) 

Considering that a typical explicit algorithm has its stability limit /^max * ^ the maximum step size 

allowed by (42) is in fact several orders of magnitude larger than allowed by any explicit algorithm. 


Computational Efficiency 

Solution of the Kalman filtering equations in second-order form is prompted by the potential gain in com- 
putational efficiency due to the beneficial nature of matrix sparsity and symmetry in the solution matrix 
of the second-order estimator equations. There is an overhead to be paid for the present second-order pro- 
cedure, that is, the additional computations introduced to minimize the control force and state estimation 
error terms on the right-hand-side of the resulting discrete equations. The following paragraphs show the 
second-order solution is most advantageous for estimator models with sparse coefficient matrices Af, D and 
K. 

Solution of the first order Kalman filter equation (6) or the second-order form (25-26, 39) may be performed 
using a time discretization as given by (19). For linear time invariant (LTI) systems, the solution matrix is 
decomposed once and subsequently upper and lower triangular system solutions are performed to compute 
the estimator state at each time step. Thus, the computations required at each time step result from 
calculation of the RHS and subsequent triangular system solutions. For the results that follow, the number 
of floating point operations are estimated for LTI systems of order 0(N). In addition, it is assumed that the 
mass, damping and stiffness matrices (M,D and K) are symmetric and banded with bandwidth aN , where 
0 < a < (0.5 - ^). 

The first-order Kalman filter equation (6) requires (AN 2 + 2Nr -b O(N)) operations at each time step. The 
discrete time second-order Kalman filter solution (25-26, 33-34) require (&a 2 N 2 +2acN 2 +3N m+ANr+0(N)) 
operations and the continuous time second-order Kalman filter (25-26, 39) require (8 aN 2 + 2aN 2 -\-5Nm + 
6iVr + (r-J- m) 2 + O(N)) operations at each time step. To examine the relative efficiency of the first-order 
and second-order forms, several cases are presented as follows. 

First, a worst case condition is examined whereby M } D and K are fully populated (c* = 0.5 — 277) and 
r =r m — N . Only for this extreme condition with large numbers of sensors and actuators relative to the 
system order, the first order Kalman filter becomes somewhat more efficient than the second-order discrete 
Kalman filter solution presented herein. 

For typical structural systems, M and K are almost always banded. In addition, the number of sensors 
and actuators is usually small compared to the system order N. If the number of actuators (m) and the 
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number of measurements (r) are proportional to the bandwidth ( r = m — ocN) t the second-order discrete 
Kalman filtering equations become computational attractive as long as a < 0.394. It should be noted that 
the larger the size of the structural systems, the smaller the bandwidth becomes, with the range of a to be 
0.05 < a < 0.15. 

Finally, for the special case of modal-based structural models, one has a — ► 0. For this case, as long as 
sensors and actuators are sufficiently smaller than the modal degrees of freedom, the present second-order 
state estimator can be substantially more efficient than the classical first-order form. This is because the 
conventional state space-based estimator must deal with a fully coupled nonsymmetric *1N x 2 N system 
whereas the present second-order estimator deals with a diagonal N x N system. A more detailed discussion 
can be found in Belvin 9 . 


Implementation and Numerical Evaluations 

The second-order discrete Kalman filtering equation derived in (25) and (26) have been implemented along 
with the stabilized form of the controller u and the filtered measurements z in such a way the estimator 
computational module can be interfaced with the partitioned control-structure interaction simulation package 
developed previously by Belvin9, Park and Belvin 17 Alvin and Park 21 . It is emphasized that the solution 
procedure of the present second-order discrete Kalman filtering equations (25) and (26) follows exactly the 
same steps as required in the solution of symmetric, sparse structural systems. It is this attribute that makes 
the present discrete filter attractive from the simulation viewpoint. For a succinct comparison between the 
present CSI simulation procedure and conventional state space-based simulation procedures, the equations 
that need to be implemented in both of the procedures are summarized below. 

Partitioned Control- Structure Interaction Equations 

The partitioned procedure for simulating the control-structure interaction problems developed in Belvin 9 and 
Park and Belvin 17 exploits the second-order diferential equation form whenever possible as shown below. 


r Structure: a) 

Sensor Output: b ) 

i Estimator: c) 

Control Force: d) 

< Estimation Error: e) 


Mq + Dq -f Kq = f 4 Bu 4 Gw 
q(0) = qo> q(°) = <io 

z = Hx 4 v 


0l f 61 I'D — Ilfql f 0 1 , [M Oil'Ll! 

1 J \ ^ + L K 0 J \qj “ (f+BuJ + [ 0 I J [ 1-2 J ^ ( ) 


M 0 

0 

5(0) = 6, P(0) = f(0) + Bu(0) 


u4F 2 M 1 Bu = F 2 (M 1 p 4 L 2 7) 4 Fi q 
74H v L 27 = z - - Bu) - 


In addition, notice that the control laws (u) and the estimation error (7) are parabolically stabilized and 
solved in a separate software module from the estimator and the structural analyzers, thus effectively ren- 
dering a computaionally efficient and accurate procedure. 
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Conventional Control-Structure Interaction Equatioons 

In contrast to the partitioned procedure summarized above, conventional control-structure interaction simu- 
lation employs a first-order differential equation form as shown below, thus requiring the solution of 2n x 2n)- 
system equations for structures and the observer. In addition, the control laws and the estimation errror are 
not stabilized, which can give rise to an accumulation of computational errors. 


where 


and 


r Structure: a) 

Sensor Output: b) 

Estimator: c) 

Control Force: d ) 


l Estimation Error: e) 


x = Ax-f Ef+ Bu + Gw 
x(0) = x 0 

z = Hx + v 

x = Ax + Ef + Bu + lor 
x(0) = 0 

u = — Fx 
7 = z — Hx 



(46) 


Numerical Experiments 

The first example is a truss beam shown in Fig. 1, consisting of 8 bays with nodes 1 and 2 fixed for 
cantilevered motions. Actuator and sensor locations, as well as their orientation, are given in Table 1. 
In the numerical experiments reported herein, we have relied on the Matlab software package 22 for the 
synthesis of both the control law gains and the discrete Kalman filter gain matrices. Figures 2, 3 and 4 
show the vertical displacement time response at node 9 for open-loop, full state feedback, and dynamically 
compensated feedback cases, respectively. In the present pepper, a full state feedback corresponds to the case 
for which the number of sensors are the same as the total system degrees of freedom whereas the dynamically 
compensated case corresopnds to a smaller number of sensors as compared with the total system degrees of 
freedom. Note the effectiveness of the dynamically compensated feedback case with four actuators and six 
sensors as indicated in Table 1 by the present second-order discrete Kalman filtering equations as compared 
with the full state feedback cases. 

Figure 5 illustrates a testbed model of an Earth-pointing satellite. For vibration control, 18 actuators and 
18 sensors are configured throughout the system; their locations are provided in Tables 2 and 3. Figures 6, 7, 
and 8 are a representative of the responses for open-loop, full state feedback, and dynamically compensated 
cases, respectively. In both examples, the estimator states are the estimated physical displacements and 
generalized momenta as previously developed, and thus the number of effective states is equal to 2Af , where 
N is the number of physical displacement variables of the second-order structural system. Therefore, the 
Kalman filter for the truss example has 108 states, and the filter for the satellite has 1164 states, a substantial 
increase over typical estimator orders for such systems. Further simulations with the present procedure should 
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shed light on the performance of dynamically compensated feedback systems for large-scale systems as they 
are computationally more feasible than heretofore possible. 

The computational overhead associated with the full state feedback vs. the use of a dynamic compensation 
scheme by the present Kalman filtering equations is reported in Table 4. It is seen that the use of the 
present second-order discrete Kalman filtering equations for constructing dynamically compensated control 
laws adds computational overhead, only an equivalent of open-loop transient analysis of symmetric sparse 
systems of order N instead of 2N x 2 N dense systems. This is evidenced in Table 4 in that the normalized 
CPU time for the dynamically compensated case (designated as K. Filter) is 63.16 whereas the total CPUs 
for the full state feedback case (FSFB) plus that of the open loop dynamic response (Transient) is 64.18. 


Summary 

The present paper has addressed the advantageous features of employing the same direct time integration al- 
gorithm for solving the structural dynamics equations and for integrating the associated continuous Kalman 
filtering equations. The time discretization of the resulting Kalman filtering equations is further facilitated 
by employing a canonical first-order form via a generalized momenta. When used in conjunction with the 
previously developed stabilized form of control laws in Park and Belvin 17 , the present procedure offers a sub- 
stantial computational advantage over the simulation methods based on a first-order form when computing 
with large (i.e. nearly full system dynamics) and sparse estimator models. 

In order to minimize the deleterious effect of numerical damping and phase distortion in the solution of the 
discrete Kalman filtering equations, the trapezoidal rule is employed. This is due to the wellknown fact that 
the trapezoidal rule conserves the system energy with minimum phase error among all the time integration 
formulas of second-order accuracy 15,16 

Computational stability of the present solution method for the filter equation has been assessed based on 
the stability analysis result of partitioned solution procedures 18 . To obtain a sharper estimate of the stable 
integration step size, a more rigorous computational stability analysis is being carried out and will be reported 
in the future. 
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Table la 


Actuator Placement for Truss Example Problem 


Actuator 

Node 

Component 

1 

2 

y 

2 

18 

y 

3 

9 

y 

4 

9 

X 


Table lb 

Sensor Placement for Truss Example Problem 


Sensor 

Type 

Node 

Component 

1 

Rate 

2 

y 

2 

Rate 

18 

y 

3 

Rate 

9 

y 

4 

Rate 

9 

X 

5 

Position 

9 

y 

6 

Position 

9 

X 



Table 2 

Actuator Placement for EPS Example Problem 


Actuator 

Node 

Component 

1 

97 

X 

2 

97 

z 

3 

96 

X 

4 

96 

z 

5 

65 

y 

6 

68 

y 

7 

59 

y 

8 

62 

y 

9 

45 

y 

10 

45 

z 

11 

70 

y 

12 

70 

z 

13 

95 

X 

14 

95 

y 

15 

95 

z 

16 

95 

4>x 

17 

95 

4>y 

18 

95 

<f>z 




Table 3 

Sensor Placement for EPS Example Problem 



3 

Rate 

96 

X 

4 

Rate 

96 

z 

5 

Rate 

65 

y 

6 

Rate 

68 

y 

7 

Rate 

59 

y 

8 

Rate 

62 

y 

9 

Rate 

45 

y 

10 

Rate 

45 

Z 

11 

Rate 

70 

y * 

12 

Rate 

70 

z 

13 

Position 

95 

X 

14 

Position 

95 

y 

15 

Position 

95 

z 

16 

Position 

95 


17 

Position 

95 

4>y 

18 

Position 

95 




Table 4 


CPU Results for ACSIS Sequential and Parallel Versions 


Model 

Problem 

Type 

Sequential 

Parallel 

54 DOF 

Transient 

4.5 

5.6 

Truss 

FSFB 

9.4 

10.2 


K. Filter 

13.0 

10.7 

582 DOF 

Transient 

98.6 

100.3 

EPS7 

FSFB 

190.2 

294.5 


K. Filter 

284.2 

321.5 





1 3 5 
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Figure 1: TRUSS BEAM PROBLEM 
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Truss Model: Open Loop Transient Response 
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Figure 2: Truss Transient Response 
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Figure 3: TRUSS FULL STATE FEEDBACK RESPONSE 
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Truss Model: Controlled Response with Kalman Filter 
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Figure 4: Truss Response with Filter 




r- 

CT> 



<o 

Oi 


Figure 5: GENERIC EARTH POINTING SATELLITE EXAMPLE 




EPS7 Model: Open Loop Transient Response 
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Figure 6: EPS TRANSIENT RESPONSE 
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EPS7 Model: Full State Feedback Response 



0.0 2.0 4.0 6.0 8.0 10. 


Time, sec 


u x at Node 45 
Uy at Node 45 
Uj at Node 45 


'igure 7: 


EPS FULL STATE FEEDBACK RESPONSE 


Model: Controlled Response w/Kalman Filter 
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Figure 8: EPS RESPONSE with KALMAN FILTER 




